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Abstract 

Background: The linear combination of variables is an attractive method in many medical analyses targeting a 
score to classify patients. In the case of ROC curves the most popular problem is to identify the linear combination 
which maximizes area under curve (AUC). This problem is complete closed when normality assumptions are met. 
With no assumption of normality search algorithm are avoided because it is accepted that we have to evaluate 
AUC n d times where n is the number of distinct observation and d is the number of variables. 

Methods: For d = 2, using particularities of AUC formula, we described an algorithm which lowered the number of 
evaluations of AUC from n 2 to n(n-1 ) + 1 . For d > 2 our proposed solution is an approximate method by considering 
equidistant points on the unit sphere in R d where we evaluate AUC. 

Results: The algorithms were applied to data from our lab to predict response of treatment by a set of molecular 
markers in cervical cancers patients. In order to evaluate the strength of our algorithms a simulation was added. 

Conclusions: In the case of no normality presented algorithms are feasible. For many variables computation time 
could be increased but acceptable. 
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Background and previous results 

In oncology one of the most used endpoint is treatment 
response. Let's denote by D the associated variable. There 
are two possible values: D = 1 if the patient responds to 
treatment and D = 0 if the patient has no response. 

Let's suppose that there are two prognostic factors and 
let's denote by X l and X 2 the random variable associated. 
Xx and X 2 could be numeric or ordinal and the patient is 
getting better or worse as the value is smaller or bigger. For 
simplicity of the talk we suppose that both are numeric. 

If X is one of X± or X 2 and c is a value from the range of 
X then the sensitivity (Se) or value "true positive" (TP) of 
variable X for c value is the probability that X > c for the 
patients which have a positive response to treatment: 

Se x (c) = TP x (c) = P(X > c\D = 1) 
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The specificity (Sp) is the probability that X < c for the 
patients which have no response "true negative" (TN): 

Sp x (c)=P(X<c\D = 0). 

A special importance has "false postive" (FP) value 
defined by 

FP x (c) = 1-Sp x (c) = P(X > c\D = 0). 

For a continuous variable X, "receiver operating 
characteristics" (ROC) curve [1] is the curve formed 
with the points 

(l-SP x (c),Se x (c)) 

that is 

(FP x (c),TP x (c)) 

for all possible values of c. 

Area under curve (AUC) "measures" the potential 
influence of the random variable on treatment response. 
AUC values are between 0.5 and 1 and if they are in the 
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proximity of 1 the variable is more important in the 
process of response prediction. 

In the case of a discrete random variable with the 
numerical values C\ < c 2 < ... < c„, the ROC curve is formed 
by joining the points 

{(0,0), (£P( Cl ), 27(a)), (FP(c 2 ), 22>(c 2 )),..., (FP(c„). TP(c„)). (1, 1)} 

For continuous variables with unknown distributions the 
simplest way to evaluate AUC is to take a random sample 
and to build the polygonal line as for discrete variables. 

The theory is similar if the signs > and < are changed 
each other in previous definitions. In practice it is chosen 
an increasing sequence c 1 <c 2 < ... <c n or a decreasing 
sequence c 1 > c 2 > ■■■ > c n so that AUC > 0.5. 

Major interest is to test equality of AUC with 0.5. 

If we have a unique random variable from all studied 
variables which has AUC > 0.5, at chosen significance level, 
than we can use this variable as prediction instrument. 

If exists multiple variables with AUC > 0.5 emerges the 
problem of multivariate prediction counting on all variables. 

Let's suppose first that we have only two random vari- 
ables. First natural variant is to choose a linear combination 
of the two variables as a global instrument of response 
prediction. 

In formal terms the problem can be stated as an 
algorithm to find a pair of real numbers [a\, a 2 ) so 
that global random variable 

Z = a\X\ + a 2 X 2 

induces a maximal AUC. 

For a clear presentation let's suppose that for the pair 
(Xi, X 2 ) there are n distinct observed values denoted by 

(x u ,x 2i ),i= 1,...,» (1) 

Also 

(n 0i ,nii),i = l,...,n (2) 

denote the number of patients that have no response, 
have response respectively for observation groups i; 

n n 
i=l i=l 

denote the whole number of patients without response, 
with response respectively. 

The ideea to solve frontal the problem without sup- 
plementary hyoptheses was generally rejected because 
at first sight the algorithms that evaluate AUC for all 
possible cases are complicated and this needs longer 
times to solve even for lower values of n and even 
with the help of computers. 

Usually this problem is solved adding supplementary 
conditions or hypoteses to variables X 1 and X 2 [2-8]. In 



[9,10] there are two comprehensive surveys. The problem 
is completly solved only when normality is supposed for 
variables X x and X 2 . As sofware we have to mention SAS 
solution of [11] for normality case. 

Present paper for a pair of variables (Xi, X 2 ) shows a 
reasonable algorithm which evaluates AUC for at most n 
(n - 1) + 1 times where n is the number of distinct 
values of the sample. For more than two variables it is 
proposed an algorithm which produces well aproximate 
solutions. 

Firstly we prove some properties of linear combinations 
of two variables which are the basis of our algorithm. Next 
paragraph introduces an approximate solution for the case 
of two and extends the algorithm to more than two 
variables. An example occured in the cancer reaserch of 
our lab is presented subsequently. The example is solved 
with programs showed in Additional file 1. For each 
program short explanations or comments are inserted. 
The paper end with a summary of a simmulation on 20 
studies with 200 observations each in order to evaluate 
the reliability of altgorithms. 

Results 

Properties of AUC evaluated for variables formed by 
linear combinations of two variables 

The algorithm from next section is based on some 
elementary proprieties derived from the calculus formula 
of AUC. 

Let's suppose that there are two real values ct\, a 2 fixed 
and we try to evaluate AUC for the linear combination 

Z = U\X\ + a 2 X 2 with the observations shown at (1) 
and (2). 

Let's denote 

Zt = cc\Xu + a 2 x 2i ,i = (4) 

the sample values of Z variable. 

From [1,9,10,12] the formula to evaluate AUC for 
random variable Z is 

1 " 

AUC = V nunojf(zi,Zj) (5) 

«o«i £=J 

where 

{1 for Zi > Zj 
0.5 for zi = Zj (6) 
0 for Zi < Zj 

with Zi,z 2 , ...,z„ sorted ascending. In practice it is chosen 
ascending or descending order of Z\, z 2 , ...,Z„ so that 
AUC > 0.5 but the results are similar. 
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Property 1 

For (ax, oc 2 ) fixed, ROC curve depends only by the 
order (increasing or decreasing) in which values Z\, 
z 2 ,...,z n are. 

Proof For fixed a.\, a 2 let's denote T{a lt a 2 ) = {zt = 
a\X\i + a 2 x 2 i\i = 1, • «} 

If r(«i, a 2 ) has m distinct elements ti < t 2 < ... < t m and 
if I t denotes the set of indexes so that the variable Z 
takes value t: I t = {i\Zi = t} then from (5) and (6) ROC 
curve depends only by the set M(ai,a 2 ) = {/ fl , ...,/(,„}. 

Property 2 

Each point located on a line through origin determines 
same ROC curves. 

Proof For a fixed pairax, a 2 , {{Xct\, Xa 2 ) |A real, A * 0} is 
the line through origin. It produces same ROC curve 
due to fact that M{Xa.\, Xa 2 ) = M(a lt a 2 ) for any A. 

Property 3 

Each pair of values i\ * i 2 with Zi, = Zi 2 determines a line 
trough origin and the points of this line generate same 
ROC curve. 

Proof Let's suppose that at least two values from the set 
T(a lt a 2 ) are equal. Let's denote i lt i 2 two indexes with 
z'i * i 2 and z^ = Zi 2 that is 

aiXi k + a 2 x 2il = aiXu 2 + a 2 x 2h 

and further 

«i (xih -xik ) + «2 {x 2il -x 2h ) = 0 (7) 

In the plane a\Qa 2 , (7) is the equation of a line that 
passes through origin. 

Property 4 

The set of points («i, a 2 ) where ROC curve has same 
value is convex. 

Proof The forth property shows that if M{a\ 1 a 2 ) =M 
j for {ai,ai)*{ol it d\ then M.(a i ,aA=M 

(ai,a 2 ) = M^a i ,a^j for any point (a i ,a^j located on 

the segment determined by {a v a 2 ) and (cc^ d^j . The 
proof comes from the observation that for any point 
(a^a^J on the segment (ai,a 2 ) and (a'^a'^J there is a 

real number A 6 [0, 1] so that a i = \a\ + (l-A)a' and 
a 2 = Xa 2 + (l-A)a 2 . We show that the order of values 

Zi,z 2 ,...,z„ remains unchanged also for ia^a]. 



Indeed, for two distinct indexes i, j with z t < Zj for both 
(ai, a 2 ) and [a i ,d^\ we compute the values of z,- and Zj for 

a , a 

1 ' T 



(Xai + (l-A)aQ#i; + (Xa 2 + {\-X)a^jx 2 i 
X(a x xii + a 2 x 2 i) + (1-A) (axu + dx^ = 
\Zi{a x , a 2 ) + (l-A)zj (V , d 2 ^ 
< Xzj(ai,a 2 ) + (1-A)z;l ot ,a 



\{a\Xij + a 2 x 2 j) + (l-Xjfa^y + oc^i^j = 
Xa\ + (1-X)a^jxy + (Xa 2 + (1-X)a^jx 2 j 
Xy + a 2 x 2 j = a 2 ). 



Further if the order of z 1 ,z 2 ,...,z n is unchanged for 
a , a ) , the ROC curves are identical. 



Algorithm to identify the linear combination of two 
variables which maximizes AUC 

We have to identify in plane (ai, a 2 ) the regions where 
AUCs are constants. From previous section we know 
that these regions are infinite triangles with the peak in 
origin. These triangles can be defined by the lines coming 

from (7). The whole number of them is C 2 n = and 
they divide the plane in maximum + 1 distinct regions. 
From the last property M{a 1} a 2 ) is constant if {a 1} a 2 ) are 
in the same region. Now we have to compute AUC for a 
point from each region and for a point from each line 
through origin that split two regions. The maximum 
number of AUC evaluations are n(n - 1) + 1. 

To finish we need a strategy to chose the points where 
AUC will be evaluated. Our proposition consists of building 
up an auxiliary line that intersects all lines (7). The intersec- 
tions with lines (7) generates maximum C„ - 1 finite seg- 
ments and two infinite segments. For the finite segments we 
have chosen the margins and the middles as points to evalu- 
ate AUC. For the infinite segments we have chosen points 
located at distance of one unit from the fixed margin. 

The authors have a program in Additional file 1 by which 
they solved the problem from above. In this program they 
have chose for auxiliary line the slope equal with 



min 



Xy X\i 
X 2 j-X 2 i 



i*j,X 2 j*X2i}-l 



(8) 



and the line passes through the point (0, 1). This slope is 
lower than all slopes derived from equations (7) so that 
the intersection points are certain. 
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Supplementary the points where we evaluate AUC can 
be chosen normalized conform to second property on 
the unity circle so that af + a 2 = 1. 

Approximate methods to identify the linear combination 
with maximal AUC 

In [5] the evaluation of («i, a 2 ) with a 2 * 0 in the 
expression U\X\ + a 2 X 2 is reduced at the identification 
of a 6 [ - 1, 1] in Xi + aX 2 and then the interval [ - 1, 
1] is divided in 201 equal segments. The maximal 
value is from the set of AUC on each segment 
extremity. Our proposition is to consider on unity 
circle all the points where AUC is evaluated. Supplemen- 
tary from symmetry we need to evaluate AUC only in 
quadrant I and IV. More exactly we evaluate AUC for 
(«i, a 2 ) with 



sinf9, a 2 ■■ 

71 n 

~ 2 '~2 ' 



cos(9for(9 

71 71 

200'" 2 + 



200' 



71 

' 2 



200- 



200' 



The precision can be improved by dividing quadrant I 
and IV in more and more regions subsequently. Practic- 
ally we divide the quadrant I and IV till the divisions are 
smaller than an apriori limit. 

This view permits easy extension when we have more 
than two prognostic factors. 

For Xi,X 2 , ...,Xf prognostic factors, with f> 2, extension 
consists in a method to highlight or to move on the 
unit sphere in space with / dimensions. Our proposal 
is to consider for a 1; a 2 , (Xf the following values: 



ai = cosd 1 

a 2 = sinf3i cos#2 

«3 = sinf?! sinf3 2 cosf9 3 

(Xf-i = sinf9i sin02... cosf?^! 
a./ = sinf^i sin02... sinf3y_! 



(9) 



with 



S.l - 7 L. 

I 2 ' 2 



TT 

200' 



TT 

2 



200 



TT 

2 



71 1 

200 }. 

200/ 



(10) 



Of course if we want to increase the precision we 
can increase the number points inside the interval 

[-!»+!]• 

The authors have a program in Additional file 1 which 
was used to solve the example from next section. 



status found at 30 days after the end of treatment. We 
have D=l or D = 0 as the patient presented complete 
remission or residual tumor at 30 days. It were 14 
patients with D = 1 and 12 patients with D = 0. 

From univariate analysis were retained: Vascular 
Endothelial Growth Factor Receptor (VEGFR) (AUC = 0.74, 
p = 0.02), dimesion of tumor (AUC = 0.73, p = 0.001) and 
age (AUC = 0.67, p = 0.06). Logistic model for multivariate 
analysis [14] did not validate any linear combination of 
these factors. 

Due to this failure we built a program associated to the 
method described in paragraph 3 (see Additional file 1). 

We started by dividing quadrant I and IV in 50 parts. 
Linear combination that maximizes the AUC for this 
division has solution: 

{0.998027, -0.0608178, 0.0156154} 
and AUC = 0.815476. 

Dividing the I-st and IV-th quadrant in 100 parts 
yields the following solution 

{0.998027, -0.0602973, 0.017518}, 
{0.998027, -0.0608178, 0.0156154}, 
{0.995562, -0.0939226, 0.00590911} 
and AUC = 0.815476. 

For 150 parts the solution is 

{0.998027, -0.0604775, 0.0168856}, 
{0.998027, -0.0608178, 0.0156154}, 
{0.996493, -0.0835127, 0.00525418} 
and AUC = 0.815476 

For 200 parts the solution is 

{0.996917, -0.0753438, 0.0218894}, 
{0.996917, -0.0756783, 0.0207032} 
and AUC = 0.821429. 

For 300 parts the solution is 

{0.997314, -0.0694128, 0.02336}, 
{0.997314, -0.0696537, 0.0226318}, 
{0.997314, -0.0698868, 0.0219012}, 
{0.997314, -0.0701124, 0.0211682}, 
{0.997314, -0.0714744, 0.0159764}, 
{0.997314, -0.0716378, 0.0152271} 
and AUC = 0.821429. 



Example 

In [13] there is an interim result of a study for several 
molecular markers in relation to response to treatment 
for cervix cancers. Endpoint was considered the patient 



As can be seen increased number of divisions for 50, 
100 and 150 does not change the maximum of AUC but 
increases the number of points where maximum AUC 
value is reached. 
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For 200 and 300 divisions the same area under the 
curve with very small increase for AUC of 0.00595238 
makes us believe that we are close to global solution. 

Figure 1 shows the ROC curves for the two linear 
combinations that give the two AUC values outlined above. 
We used firstly the score: 

0.998027 x VEGFR-0.0608178 x dimension of tumor 
+ 0.0156154 x age 

resulting from algorithm with 150 divisions then the 
score 

0.996917 x VEGFR-0. 0753438 x dimension of tumor 
+ 0.0218894 x age 

resulting from algorithm with 300 divisions. 

Note that both scores have the values of p highly 
significant, and we propose the solution that has 
higher AUC. 

Although computer times were acceptable (between 
37 seconds to 50 divisions and 1 hour and 13 minutes to 
300 divisions) we do not believe that would be necessary 
to go further with the number of divisions and we believe 
that a good solution could be 

0.996917 x VEGFR-0.0753438 x dimension of tumor 
+ 0.0218894 x age. 

Furthermore criteria of classification from ROC curve 
analysis [9] tells us for this choice that patients with 



score higher than 1.425782 are patients from whom we 
expect a better result (Se = 0.71, Sp = 0.92). 

Simulation 

As previous example has a small number of observations 
we have made a simulation for 20 studies with 200 
observations each with three prognostic factors. For the 
first factor, cases were selected from a pseudonormal 
variable with mean 1 and standard deviation of 3 and 
controls from a pseudonormal variable with mean 3 and 
standard deviation 3.5. The second and third prognostic 
factor, also come from a pseudonormal variable with 
standard deviation of 3 and 3.5 respectively for cases 
and controls and with averages of 4 and 6 for controls 
respectively 6 and 6.5 for cases. 

The simulation was made on a Lenovo computer with 
operating system Windows 7 Ultimate on 64-bit with an 
i7 processor at 1.37 Gz in parallel with the current work 
of the author i.e. text editing, Internet browsing and 
reading emails. The result of the simulation for the 
algorithm presented before for 50, 100 and 200 segments 
are shown in Table 1. It is noted that the jump from 
50 to 100 segments produces a change in AUC only 
to the third decimal place (the maximum value of 
0.0012 to simulation 15). Jump from 100 segments to 
200 segments changes AUC only at the fourth decimal 
place (the maximum value of 0.0008 to simulation 14). 
We believe that in practice there is no need to move 
beyond 200 divisions only for outstanding situations. 



i.o - 



5.8 



C.6 



0.4 - 



0.2 - 
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02 



0.4 



C.6 



C.S 



1.0 
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Figure 1 ROC curves for score 0.998027 xVEGFR -0.06081 78 x dimension of tumor + 0.01 561 54 x age with AUC = 0.81 5476 and 

p = 0.000093 for 150 divisions (continuos line) and for score 0.99691 7 xVEGFR- 0.0753438 x dimension of tumor + 0.021 8894 x age with 

AUC = 0.821429 and p = 0.000056 for 300 divisions (dashed line). 
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Table 1 Results of 20 simulations with 200 observations 



Crt.Nb. 


Time 


AUC50 


Time 


AUC100 


Time 


AUC200 


AUC100 - AUC50 


AUC200- AUC100 


1 


1314s (OH 21 M 54 s) 


0.7091 


5033 s (1 H 23M 53 s) 


0.7098 


20354 s (5H 39M 14 s) 


0.7098 


0.0007 


0.0000 


2 


1283s (OH 21 M 23 s) 


0.6589 


5154 s (1H 25M 54 s) 


0.6589 


31636 s (8H 47M 16 s) 


0.6589 


0.0000 


0.0000 


3 


1501s (OH 25M 1 s) 


0.6406 


5842 s (1H 37M 22 s) 


0.6412 


23352 s (6H 29M 12 s) 


0.6412 


0.0006 


0.0000 


4 


1173s (OH 19M 33 s) 


0.6862 


4681 s (1H 18M 1 s) 


0.6862 


25012 s (6H 56M 52 s) 


0.6867 


0.0000 


0.0005 


5 


1277s (OH 21 M 17 s) 


0.6629 


10790 s (2H 59M 50s) 


0.6633 


12321 s (3H 25M 21 s) 


0.6638 


0.0004 


0.0005 


6 


1353s (OH 22M 33 s) 


0.6715 


4574 s (1H 16M 14 s) 


0.6717 


15292 s (4H 14M 52 s) 


0.6726 


0.0002 


0.0009 


/ 


1342s (OH 22M 22 s) 


0.6761 


5132 s (1H 25M 32 s) 


0.6772 


18625 s (5H 10M 25 s) 


0.6773 


0.001 1 


0.0001 


8 


1297s (OH 21 M 37 s) 


0.6944 


6988 s (1H 56M 28 s) 


0.6953 


18813 s (5H 13M 33 s) 


0.6954 


0.0009 


0.0001 


9 


1070s (OH 17M 50s) 


0.6988 


5399 s (1H 29M 59 s) 


0.6990 


19498 s (5H 24M 58 s) 


0.6994 


0.0002 


0.0004 


10 


536 s (OH 8M 56 s) 


0.6638 


3022 s (OH 50M 22 s) 


0.6640 


18556 s (5H 9M 16 s) 


0.6646 


0.0002 


0.0006 


11 


1329s (OH 22M 9 s) 


0.6900 


4766 s (1H 19M 26 s) 


0.6902 


20419 s (5H 40M 19 s) 


0.6906 


0.0002 


0.0004 


12 


1288s (OH 21 M 28 s) 


0.6946 


5086 s (1H 24M 46 s) 


0.6948 


20573 s (5H 42M 53 s) 


0.6948 


0.0002 


0.0000 


13 


637 s (OH 10M 37 s) 


0.6873 


2454 s (OH 40M 54 s) 


0.6875 


21271 s (5H 54M 31 s) 


0.6875 


0.0002 


0.0000 


14 


513 s (OH 8M 33 s) 


0.7031 


2025s (OH 33M 45 s) 


0.7032 


20139 s (5H 35M 39 s) 


0.7040 


0.0001 


0.0008 


15 


952 s (OH 15M 52 s) 


0.7200 


2082s (OH 34M 42 s) 


0.7202 


21224 s (5H 53M 44 s) 


0.7204 


0.0002 


0.0002 


16 


1176s (OH 19M 36 s) 


0.7401 


4923 s (1 H 22M 3 s) 


0.7413 


27836 s (7H 43M 56 s) 


0.7413 


0.0012 


0.0000 


17 


796 s (OH 1 3M 1 6 s) 


0.7398 


4332 s (1H 12M 12 s) 


0.7399 


18213 s (5H 3M 33 s) 


0.7405 


0.0001 


0.0006 


18 


1296s (OH 21 M 36 s) 


0.6635 


2534 s (OH 42M 14 s) 


0.6638 


20165 s (5H 36M 5 s) 


0.6644 


0.0003 


0.0006 


19 


797 s (OH 13M 17 s) 


0.7041 


3407 s (OH 56M 47 s) 


0.7045 


20313 s (5H 38M 33 s) 


0.7045 


0.0004 


0.0000 


20 


1420s (OH 23M 40s) 


0.6825 


5532 s (1H 32M 12 s) 


0.6826 


15051 s (4H 10M 51 s) 


0.6826 


0.0001 


0.0000 


Average 


1117s (OH 18M 37 s) 




4687 s (1 H 1 8M 7 s) 




20433 s (5H 40M 33 s) 









AUC50, AU CI 00, AUC200 denotes the approximation of AUC by dividing interval [-§, + §] in 50, 100, 200 equal parts. 



Average execution time for 50, 100 and 200 segments 
was 18 minutes, 1 hour and 18 minutes, 5 hours and 
40 minutes which is an acceptable time for a practical 
problem. 

Discussions and conclusions 

Multivariate analysis is used largely in any medical 
paper. However testing the hypotheses in modeling is 
not a very simple task and this is the reason for trying a 
lot of potential models and choose the model best suited 
to observations. The papers of [4,6,15-19] prove that 
there is a large basis to use linear combinations of 
variables in ROC analysis. If we do not have solid 
condition to apply for example one of the cited 
models, the method from our paper produces always 
a score for which we have maximal AUC or an 
approximate. 

On the other hand in a classical model of regression it 
is known that the numerical methods used to identify 
the model parameters not always provide a global 
maximum and depends heavily on the initial values of the 
algorithm. The solution presented we believe could be 
used there as a baseline for these algorithms. 

The main advantage of presented algorithms is that 
it always provides a solution. However for many prognostic 



factors and observations, time of the calculation could be a 
problem. 

Certainly, approximate method is more appropriate in 
this last case despite the fact that it does not guarantee a 
global solution. However it is guaranteed to yield a solution 
with AUC higher than each variable taken separately. 

Our algorithm can be used in any medical paper as an 
alternate method for multivariate analysis. 

The presented algorithms have major advantage to pro- 
vide always a solution with no supplementary constraints. 

For many variables computation time is high but not 
high enough as not to accept this cost. 

Additional file 



Additional file 1: AUC evaluation, Maximum AUC evaluation for a pair 
of variables, Maximum AUC evaluation for more than two variables. 
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